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Methane and carbon dioxide adsorption on edge-functionalized graphene: A 
comparative DFT study 
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With a view towards optimizing gas storage and separation in crystalline and disordered nanoporous carbon- 
based materials, we use ah initio density functional theory calculations to explore the effect of chemical 
functionalization on gas binding to exposed edges within model carbon nanostructures. We test the geometry, 
energetics, and charge distribution of in-plane and out-of-plane binding of CO2 and CH4 to model zigzag 
graphene nanoribbons edge-functionalized with COOH, OH, NH2, H2PO3, NO2, and CH3. Although different 
choices for the exchange-correlation functional lead to a spread of values for the binding energy, trends across 
the functional groups are largely preserved for each choice, as are the final orientations of the adsorbed gas 
molecules. We find binding of CO2 to exceed that of CH4 by roughly a factor of two. However, the two 
gases follow very similar trends with changes in the attached functional group, despite different molecular 
symmetries. Our results indicate that the presence of NH2, H2PO3, NO2, and COOH functional groups can 
significantly enhance gas binding with respect to a hydrogen-passivated edge, making the edges potentially 
viable binding sites in materials with high concentrations of edge carbons. To first order, in-plane binding 
strength correlates with the larger permanent and induced dipole moments on these groups. Implications for 
tailoring carbon structures for increased gas uptake and improved CO2/CH4 selectivity are discussed. 
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I. BACKGROUND 



Nanoporous carbon-based materials are particularly 
attractive as gas storage media, due to their high in- 
trinsic adsorptive capacity, low cost and weight, and rel- 
ative ease of processing. Accordingly, both crystalline 
and disordered nanoporous carbons have been the sub- 
ject of active investigation for vehicular storage of nat- 
ural gas and hydrogen, as well as framework materials 
for CO2 capturcii"— Among the adjustable processing 
parameters is the ability to tune the internal chemistry 
via insertion of desired functional groups to enhance the 
interaction between the substrate and gas adsorbate or 
increase the active surface area of the material, thereby 
increasing uptake. If the substrate-adsorbate interaction 
can be enhanced selectively, the host material may also 
be used for industrial processes such as natural gas pu- 
rification, where efficient separation of CO2 contaminants 
from CH4 is necessary. ^^"^^ In disordered porous materi- 
als, such as activated carbons and carbon aerogels^^ii^ 
this is often achieved by thermal activation, pyrolysis, or 
chemical/electrochemical treatment. Depending on the 
specific nature of the chemical process employed, this can 
result in a predominance of certain classes of functional 
groups ji2i2£ The intrinsic disorder gives rise to significant 
local deviations from the dominant planar sp^ geome- 
try, including reasonable concentrations of reactive edges 
that act as ready binding sites for additional functional 
groups fi^^"— Efforts towards increasing gas uptake and 
selective separation have often focused on understanding 
the influence of activation on the resulting pore size dis- 
tribution and available surface area. Generally speaking. 



less attention has been paid to the specific interactions of 
the functional groups with the adsorbate. This is in part 
due to complexities in characterizing which groups are 
present, as well as in separating the effects of increased 
surface area and topological defects from those of chem- 
ical functionalization.^^ 

Crystalline porous framework materials, such as metal- 
organic frameworks (MOFs) and covalent-organic frame- 
works (COFs), offer additional tight control over the local 
chemistry within the framework itself. This is achieved 
by direct substitution or functionalization of molecular 
building blocks, which commonly feature segments of in- 
terconnected aromatic rings. ■^^^1^"— 

For both crystalline and disordered nanoporous carbon 
substrates, it is highly desirable to devise specific, well- 
founded strategies for modifying gas uptake via chemical 
functionalization or substitution. However, doing so re- 
liably requires a robust understanding of how changes 
in the identity and local geometry of functional groups 
attached to edges of sp^ carbon can affect the interac- 
tion with adsorbed gas molecules. Accordingly, we have 
used ab initio density functional theory to systemati- 
cally study the influence of certain commonly occurring 
edge functional groups2ii2^i2£ on methane and carbon 
dioxide adsorption within model carbon structures. The 
substrates are modeled as edge-functionalized graphene 
nanoribbons, which represent the simplest extended sp^ 
carbon system with exposed edges, and as such are suit- 
able for exploring edge functionalization in locally planar 
carbon structures 1^1^ Our results indicate that it is pos- 
sible to activate exposed edges to act as additional gas 
binding sites, thereby increasing the active surface area 
and total gas uptake of the material. Given current in- 



terest towards more accurate theoretical treatments of 
the weak dispersion forces that dominate gas binding on 
these substrates, we also provide a careful comparison 
of values obtained using recently developed exchange- 
correlation functionals. We conclude that in most cases, 
trends are well preserved across competing methods. 



II. COMPUTATIONAL DETAILS 

We use density functional theory (DFT) to determine 
the geometry and energetics of CH4 and CO2 adsorbed 
on zigzag graphene nanoribbons (ZGNRs) with attached 
chemical functional groups. We test three electron- 
donating groups (OH, NH2, CH3) and three electron- 
withdrawing groups (NO2, COOH, H2PO3). Lowdin 
population (LP) analysis is employed to analyze the 
charge redistribution resulting from adsorption. We have 
restricted ourselves to considering chemical modification 
of only the zigzag edges of graphene nanoribbons; this 
choice is motivated by the results of earlier studies show- 
ing that zigzag edges are more reactive than armchair 
edges, and are therefore more likely binding sides for edge 
functional groups. •^^'•^■* The reactivity of the zigzag edges 
was verified in our own test calculations, which indicated 
that unpassivated ZGNR edges can react spontaneously 
with nearby gas adsorbates. Similar behavior was re- 
ported in Ref. i35i. 

As physisorption dominates the interactions between 
the adsorbate and substrate, it is important to select 
a method of treating van der Waals (vdW) dispersion 
forces within a DFT framework that will qualitatively 
preserve the trends we are attempting to describe. This 
remains an area of active research;^ For methane ad- 
sorption, we compare three choices for the exchange- 
correlation functional: 1) the uncorrected Perdew-Zunger 
local density approximation (DFT-LDA)j^ 2) the local 
atomic potential formulation (DFT-|-LAP)i2^ and 3) the 
van der Waals density functional (vdW-DF),^>^iia The 
DFT-I-LAP and vdW-DF functionals are specifically de- 
signed to reproduce accurate binding energies for ph- 
ysisorbed systems by taking into account van der Waals 
dispersion forces. 

The DFT-I-LAP functional is based on the augmenta- 
tion of an existent pseudopotential by an additional chan- 
nel, expressed in an atom-centered local form,^^ The cor- 
rection is constant for short-range interactions but falls 
off exponentially for long-range interactions. The three 
independent fitting parameters in DFT-I-LAP are de- 
rived from binding curves of small dimer systems, as ob- 
tained from coupled-cluster calculations including single, 
double, and perturbative triple excitations (CCSD(T)). 
Here, we apply the correction to pseudopotentials gener- 
ated using the revPBE functional4i It should be noted 
that DFT-I-LAP values are not given for the H2PO3- 
functionalized ZGNR due to the unavailability of a pa- 
rameterized local correction for phosphorus. 

The vdW-DF functional mixes short-range LDA cor- 



relation with a fully nonlocal, ab initio approximation 
of long-range correlation, which contains no empirical 
parametersi^ In principle, the correction may be im- 
plemented self-consistently. However, numerous calcu- 
lations of weak-binding molecular systems have shown 
that non-self consistent perturbation of a properly cho- 
sen exchange functional is usually sufficient, resulting 
in substantial computational savings. ■^^''*^. Accordingly, 
we adopt the non-self consistent formulation of vdW-DF 
here. As with DFT-I-LAP, we have used the revPBE 
functional^ to compute the reference exchange term. 

All DFT calculations were performed using the 
Quantum-ESPRESSO code4^ For the DFT-LDA 
and DFT-I-LAP calculations, we used ultrasoft 
pseudopotentials^!^ with plane-wave cutoffs of 30 Ry 
and 360 Ry for the wavefunctions and charge density, 
respectively. For the vdW-DF calculations, norm- 
conserving pseudopotentials with a plane-wave cutoff 
of 80 Ry were used. For the periodic ribbon calcula- 
tions, a vacuum spacing of 20 A was introduced along 
non-repeating directions, and a 6 x 1 x 1 Monkhorst- 
Pack k-point mesh was used. Marzari-Vanderbilt cold 
smearing^ with a width of 0.007 Ry was introduced 
for improved convergence and for savings in k-point 
density. In most cases, calculations were performed with 
unrestricted spin; however, for the vdW-DF calculations, 
spin-restricted DFT was used, since the methodology is 
not rigorously defined for the spin-unrestricted case. In 
the spin-restricted calculations, all functionalized ZGNRs 
demonstrate semiconducting or half-semiconducting be- 
havior, in agreement with previous theoretical work on 
narrow ribbonsi^i^i^ In addition, we observe that the 
ground-state ZGNR spin configuration consists of two 
sets of ferromagnetically aligned edge states which are 
oppositely oriented at opposite ribbon edges, as has 
been reported elsewhere . ^^'^^ 

For comparison with the DFT adsorption numbers on 
the functionalized ZGNRs, second-order M0ller-Plesset 
perturbation theory (MP2) adsorption calculations were 
performed for methane on similarly functionalized ben- 
zene rings. MP 2 calculations are known to capture cor- 
relation effects, including van der Waals dispersion. We 
used the Dunning aug-cc-pVDZ basis set^ to describe 
all atoms, with the counterpoise correction method^ 
employed to account for basis set superposition error. 
MP2 calculations were performed using the Gaussian09 
software. ^^ Isolated benzene rings were substituted for 
the ZGNRs in the MP 2 calculations in order to make 
the calculations computationally accessible. Although 
the lower coordination of the benzene rings results in 
weaker overall adsorption and precludes direct quanti- 
tative comparison against our ZGNR values, the MP2 
binding trends across the functional groups are useful as 
an additional benchmark. They also aid in forming a 
qualitative assessment of the relative contribution of the 
functional group versus the ribbon in determining ad- 
sorption behavior. 

The model ZGNRs used in our calculations are six car- 



bon atoms wide with a single functional group attached 
to a carbon along one zigzag edge. Since we are inter- 
ested in dilute functionalization, we include in the super- 
cell three additional zigzag-edge carbons adjacent to the 
functional group, resulting in a total of 24 carbon atoms 
(see Fig. [1]). These additional carbons are terminated 
with hydrogen, as are all dangling bonds on the ribbon 
edge opposite to the functional group. As a reference, gas 
adsorption is also studied on a fully hydrogen-passivated 
ZGNR. The benzene molecules studied using MP2 also 
have a single attached functional group, with hydrogens 
attached to the remaining carbons in the ring. The ad- 
sorption energy is defined as i^ads = — ^-g+s + Eg + E^, 
where Eg+s is the energy of the adsorbed system. Eg is 
the energy of an isolated gas molecule (CH4 or CO2), and 
Es is the energy of the substrate alone. 

In order to better span the configurational phase space 
of possible CH4 and CO2 geometries, we examine two 
different initial adsorbate binding sites for each complex, 
from which the system is allowed to relax without con- 
straint. These are shown in Fig.[lja). In the first, which 
we refer to as "in-plane" (IP) adsorption, the CII4 or 
CO2 molecule is initially placed adjacent to the func- 
tional group, with the carbon atom of the gas molecule 
lying in the basal plane. In the second, referred to as 
"out-of-plane" (OP) adsorption, the gas molecule is ini- 
tially placed directly above the functional group, perpen- 
dicular to the basal plane. These two configurations were 
selected based on symmetry considerations so as to max- 
imize (OP) or minimize (IP) interaction with the ZGNR 
or benzene simultaneously with the functional group. 
Notably, these same sites were found as the two local 
extrema in a recent comprehensive study of the poten- 
tial energy surface of water binding to benzene, further 
motivating their choice herc i^^'^^ 

Three possible high-symmetry orientations exist for 
bound methane, depicted schematically in Fig. [TJb) . In 
the first, which we label "tripod toward" (TT), one of 
the methane C-H bonds along the C3 molecular axis 
points toward the functional group. In the second, "tri- 
pod away" (TA), one of the C-H bonds along the C3 
axis points away from the functional group. In the 
third, "straddle" (S), the C2 axis aligns with the func- 
tional group. Similarly, CO2 can be found in two high- 
symmetry orientations, shown in Figlljc). The first, 
labeled "axis toward" (AT), has the Coo axis of CO2 
aligned with the functional group. The second, labeled 
"axis perpendicular" (AP), instead has the C2 axis of the 
CO2 aligned with the functional group. Each of these 
was tried as the starting geometry, although we report 
only the lowest-energy final configurations here. When 
recording the final configurations for OP binding, we also 
indicate the approximate final location of the adsorbed 
gas molecule with respect to the basal plane, either above 
one of the functional group atoms (G) or else above the 
carbons on the ribbon edge (E). 

Given the difficulty of controlled synthesis of edge- 
functionalized graphene, we cannot directly benchmark 
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FIG. 1. Schematic representation of the tested adsorption 
geometries listed in Tables U and |II] and in the text. Here, 
R represents the functional-group adsorption site. Panels il- 
lustrate (a) in-plane (IP) and out-of-plane (OP) adsorption; 
(b) CH4 adsorption oriented with the tripod towards (TT) 
the functional group; (c) CH4 adsorption oriented with the 
tripod away (TA) from the functional group; (d) CH4 ad- 
sorption oriented such that two C-H bonds straddle (S) the 
functional group; (e) CO2 adsorption oriented with the C=0 
axis toward (AT) the functional group; and (f) CO2 adsorp- 
tion oriented with the C=0 axis perpendicular (AP) to the 
functional group. 



our nanoribbon results against experiments. Accord- 
ingly, we have also calculated OP binding of CH4 and 
CO2 to clean graphene sheets in a 4 x 4 unit cell, which 
allows for a more straightforward comparison. Our val- 
ues of 12.1 kJ/mol for DFT-LDA and 11.1 kJ/niol for 
DFT+LAP compare favorably with the experimental 
value of 12.2 kJ/mol for methane on graphite,— and 
match well with previous DFT-LDA calculations.— The 
effect on gas binding of adding layers to single-layer 
graphene was found to be negligibly small, both in our 
test calculations and in previous studies j^>^ It is worth 
noting that DFT-LDA reproduces experimental methane 
binding on graphite quite well without correction, in 
agreement with previous analyses of similar systems in 
the literature..^"-^ Our OP vdW-DF value of 16.9 kJ/mol 
is slightly larger than the value of 15.4 kJ/mol in Ref . ISa. 
although both numbers are higher than the experimen- 
tal value, in keeping with reports of overestimation of 
vdW-DF binding on tt systems.— We also benchmarked 
our MP2 calculations against available CCSD(T) calcu- 
lations for methane binding above the ring centers of 



benzene and phenol;^ our calculated binding energies 
of 6.4 kJ/mol and 6.5 kJ/mol for benzene and phenol, 
respectively, are similar to but slightly higher than the 
CCSD(T) values of 5.9 and 6.2 in the complete basis 
set limit. This is consistent with similar analyses in the 
literature.-^SiSi 



III. METHANE ADSORPTION 

Adsorption energies for IP and OP methane binding to 
the functionalized ZGNRs are given in Fig. [5] and Table |T] 
for three different DFT exchange-correlation functionals: 
DFT-LDA, DFT+LAP, and vdW-DF, as weU as MP2 
results for methane on similarly functionalized benzene 
rings. It is immediately obvious from Fig. [5] and Table |T] 
that the different techniques can lead to a wide range of 
values for i^ads- This is particularly true for IP adsorp- 
tion, where DFT-LDA overbinds by about a factor of two 
with respect to the van der Waals-corrected DFT meth- 
ods. For OP adsorption, the energies for the ZGNRs are 
much closer, usually falling within 30% of one another. 
As expected, addition of van der Waals corrections to 
DFT decreases E'ads for IP binding, compensating the 
overbinding within DFT-LDA. However, in the case of 
vdW-DF, these same corrections increase i?ads relative to 
DFT-LDA for OP binding. Note that vdW-DF has been 
shown to overbind tt systems, ^"^ which may explain why it 
gives larger values for OP i^ads compared with the other 
methods. For both geometries, the MP2-derived i?ads is 
significantly weaker, reflecting the lower coordination of 
the benzene system. 

Despite the quantitative differences, the trends in -Eads 
for IP adsorption across the various functional groups are 
very similar for the three DFT functionals, as well as the 
MP2 calculations (Fig. [2]). Interestingly, of the DFT re- 
sults, the trend for uncorrected DFT-LDA is actually the 
closest match to the trend for MP2, to within a multi- 
plicative constant (see Fig. [2Ic)). This finding is con- 
sistent with previous results on physisorbed systems, for 
which DFT-LDA trends compared favorably with com- 
peting exchange and correlation functionals^^-S^ and to 
Quantum Monte Carlo calculations of weak molecular 
binding to benzene<^i^ Similarly, DFT-LDA has also 
been shown to reproduce the correct preferred orienta- 
tions for methane crystals?^ and to provide an accurate 
description of electron density in weak binding of small 
molecules to carbon systems <^i^ The greatest deviations 
in the computed trends are found in the vdW-DF descrip- 
tion of IP binding to the two largest functional groups, 
COOH and H2PO3 (Fig.IHc)). Binding to these appears 
relatively weaker than one would expect by comparison 
with the other DFT methods or with the MP2 calcula- 
tions. It is unclear if this is an effect of the size of these 
groups or if it is related to the specific description of the 
methane binding orientation, which involves simultane- 
ous interaction with an O-H bond and an oxygen lone 
pair. 



For methods other than DFT-LDA, OP adsorption 
is preferred over IP adsorption. However, the differ- 
ent methods exhibit less variation across the functional 
groups for OP adsorption than for IP adsorption, mak- 
ing it more difficult to define specific trends. Moreover, 
a few of the groups did not have stable OP binding sites 
in certain descriptions, but rather collapsed into the IP 
configuration or else migrated to a faraway site above the 
ribbon center. Note that the existence or absence of sta- 
ble OP binding sites is method dependent; for instance, 
vdW-DF and MP2 identify a stable OP minimum for 
OH, whereas DFT-LDA and DFT+LAP do not. This 
suggests a qualitative dependence of the potential en- 
ergy surface for OP adsorption on the chosen electronic 
structure description, which contrasts with what we ob- 
serve for IP adsorption. An additional feature unique to 
OP adsorption is the near constancy of the MP2 results 
across the functional groups. This most likely reflects the 
differences between small benzene rings and nanoribbons. 
Nanoribbons can accommodate additional charge reorga- 
nization from functionalization through their extended tt 
manifold. 

Table U also lists the final orientations of the IP- and 
OP-adsorbed methane for each of the tested methods 
and functionalizations. These orientations are depicted 
in Fig. [3] and Fig. H for DFT-hLAP. For both IP and OP 
adsorption, the orientations are very consistent across the 
DFT methods, with tripod-toward (TT) by far the most 
common. The MP2 results are generally compatible with 
the DFT orientations, although for three of the groups 
(OH, NO2, and unfunctionalized benzene), the straddle 
(S) orientation is preferred over the tripod-toward (TT) 
orientation for IP binding. The availability of nearby 
carbon rings for direct interaction with methane C-H 
bonds leads to a TT orientation for methane on graphene, 
whereas a TA orientation is preferred on benzene. 

The specific location of the adsorbed methane for 
OP binding is also listed in Table |T] and depicted for 
DFT+LAP in Fig. H The resuhs show that adsorp- 
tion usually occurs near the zigzag edge carbon to which 
the functional group is anchored, rather than above the 
atoms of the functional group itself. The only exceptions 
to this are the unfunctionalized ZGNR, where DFT-LDA 
and DFT+LAP show stable binding above the hydrogen 
edge atoms; and the DFT-LDA result for H2PO3, which 
shows stable binding approximately above the hydrogen 
atom in the 0-H bond (this is the likely origin of the sig- 
nificant enhancement in i?ads with respect to the other 
methods, where a different minimum was found). This 
suggests that the functional group changes the electronic 
character of the carbon atom that couples it to the ex- 
tended TT system. In addition, the OP geometry features 
the C-H bond that is farthest from the functional group 
aligned so as to point toward a ring center (see Fig. \^. 
This echoes the stable binding configuration for methane 
on graphene, where each C-H bond orients toward a ring 
center (Fig. |3Kh)). 

Comparison of the binding distances for the various 
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FIG. 2. Trend in calculated adsorption energies -Bads for (a) IP and (b) OP binding of methane on ZGNRs across various 
functional groups using DFT-LDA (stars), DFT+LAP (circles), and vdW-DF (crosses) exchange-correlation functionals. 'None' 
refers to the fully hydrogen-passivated ZGNR. Adsorption energies are also shown for MP2 calculations (triangles) of methane 
binding to similarly functionalized benzene rings. Panel (c) shows the correlation of the DFT results on the ZGNRs with the 
MP2 results on benzene for IP binding. 
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FIG. 3. Final adsorption geometries for IP adsorption of 
methane on (a) an unfunctionalized ZGNR; as well as ZGNRs 
functionahzed with (b) OH, (c) NH2, (d) COOH, (e) CH3, (f) 
NO2, and (g) H2PO3. All shown configurations represent the 
DFT-I-LAP resuh, except the H2P03-functionalized ZGNR, 
which represents the vdW-DF result. Additional interactions 
of the methane with nearby passivation hydrogens are evident 
in (a), (b), (c), and (f). Color scheme: H (cyan), C (yellow), 
O (red), N (green), P (orange). The methane carbon is shown 
in gray. 
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FIG. 4. Final adsorption geometries for OP adsorption of 
methane on (a) an unfunctionalized ZGNR; as well as ZGNRs 
functionahzed with (b) OH, (c) NH2, (d) COOH, (e) CH3, 
(f) NO2, and (g) H2PO3; and (h) on clean graphene. All 
shown configurations represent the DFT-I-LAP result, except 
the H2PO3- and OH-functionalized ZGNRs, which represent 
the vdW-DF result. Color scheme: H (cyan), C (yellow), O 
(red), N (green), P (orange). The methane carbon is shown 
in gray. 



TABLE I. Final relaxed geometries and adsorption energies -Eads of CH4 adsorbed on functionalized and unfunctionalized 
ZGNRs. The final orientation and position of the CH4 molecule is given alongside the adsorption energy in parentheses (see 
text for description of acronyms). For IP adsorption, dmin represents the distance from the CH4 carbon atom to the nearest 
functional-group atom, whose identity is given in parentheses. For OP adsorption, d^ is the distance of the adsorbate from the 
ribbon basal plane. "Unstable" refers to configurations for which no stable minimum was found. 



Group 


Site 


Quantity 


DFT-LDA 


DFT+LAP 


vdW-DF 


MP2 (benzene) 


None (H) 


IP 


-Eads (kJ/mol) 


8.4 (TT) 


4.5 (TT) 


6.8 (TT) 


3.6 (TT) 






dmin (A) 


2.59 (H) 


2.85 (H) 


3.03 (H) 


3.05 (H) 




OP 


-Eads (kJ/mol) 


8.9 (TT,G) 


5.3 (TT,G) 


Unstable 


Unstable 






d. (A) 


2.67 


3.05 


Unstable 


Unstable 


OH 


IP 


Eads (kJ/mol) 


15.9 (TT) 


8.9 (TT) 


8.9 (TT) 


5.5 (S) 






dmin (A) 


2.23 (H) 


2.51 (H) 


2.89 (H) 


2.58 (H) 




OP 


E;ads (kJ/mol) 


Unstable 


Unstable 


11.2 (TT,E) 


7.4 (TT,E) 






d. (A) 


Unstable 


Unstable 


3.56 


3.56 


NH2 


IP 


Eads (kJ/mol) 


14.2 (S) 


9.8 (TT) 


8.7 (TT) 


4.0 (TT) 






dmin (A) 


2.47 (H) 


2.64 (H) 


3.11 (H) 


2.93 (H) 




OP 


-Eads (kJ/mol) 


12.4 (TT,E) 


12.8 (TT,E) 


14.9 (TT,E) 


7.6 (S,E) 






d. (A) 


3.18 


3.51 


3.72 


3.53 


COOH 


IP 


Eads (kJ/mol) 


18.2 (TT) 


9.9 (TT) 


7.2 (TT) 


5.4 (TT) 






dmin (A) 


2.08 (H) 


2.36 (H) 


2.70 (H) 


2.58 (H) 




OP 


E;ads (kJ/mol) 


13.1 (TT,E) 


11.4 (TT,E) 


15.5 (TT,E) 


7.8 (TT,E) 






d. (A) 


3.20 


3.43 


3.72 


3.46 


CH3 


IP 


E^ads (kJ/mol) 


5.9 (TT) 


2.8 (TT) 


5.5 (TT) 


2.3 (TT) 






dmin (A) 


2.95 (H) 


2.99 (H) 


3.47 (H) 


3.47 (H) 




OP 


E^ads (kJ/mol) 


11.8 (TT,E) 


9.9 (TT,E) 


Unstable 


7.8 (TT,E) 






d. (A) 


3.36 


3.56 


Unstable 


3.65 


NO2 


IP 


Eads (kJ/mol) 


12.7 (S) 


9.1 (S) 


10.6 (S) 


4.5 (S) 






dmin (A) 


3.28 (O) 


3.57 (O) 


3.92 (0) 


3.97 (O) 




OP 


Sads (kJ/mol) 


14.0 (TT,E) 


11.4 (TT,E) 


15.4 (TT,E) 


7.4 (S,E) 






d. (A) 


3.22 


3.45 


3.72 


3.55 


H2PO3 


IP 


Eads (kJ/mol) 


22.0 (TT) 


N/A 


9.4 (TT) 


7.2 (TT) 






dmin (A) 


2.04 (H) 


N/A 


2.64 (H) 


2.55 (H) 




OP 


-Eads (kJ/mol) 


22.2 (TT,G) 


N/A 


15.9 (TT,E) 


8.4 (TT,E) 






d. (A) 


3.18 


N/A 


3.64 


4.51 


Graphene/ 


OP 


-Eads (kJ/mol) 


12.1 (TT) 


11.1 (TT) 


16.9 (TT) 


6.4 (TA) 


Benzene 




d. (A) 


3.17 


3.33 


3.58 


3.73 


(ring center) 















methods in Table H] reveals that distances for DFT- 
LDA are generally shortest, followed by DFT-I-LAP, with 
vdW-DF the longest. The experimental binding dis- 
tance of methane to graphite, 3.45A,^'^ falls between 
the DFT-I-LAP (3.33A) and vdW-DF (3.58A) results for 
graphene. It is worth noting that based on Table HI bind- 
ing distances alone cannot be considered a reliable pre- 
dictor of -Eads- 

From the final IP relaxed geometries in Fig. [31 it can be 
seen that methane aligns primarily along the 0-H bond 
for the OH, COOH, and H2PO3 functionalizations, along 
the N-H bond for NH2, and along the N-0 bond for NO2. 
The impact of these polar bonds on methane binding 
geometry suggests significant interaction with the per- 



manent dipole of the functionalized substrate; the addi- 
tional importance of this quantity in predicting adsorp- 
tion strength is explored below. In addition, the smaller 
polar functional groups feature significant secondary in- 
teractions with one of the neighboring C-H bonds along 
the zigzag edge. This is the case for OH, NH2, and NO2. 
This enables the methane to orient with its tripod to- 
ward (TT) the polar bond in the functional group while 
simultaneously interacting with the side C-H bond in a 
straddle (S) orientation (these are reversed for NO2, with 
S and TT interactions for the N-0 and C-H bonds, re- 
spectively). 

The larger COOH and H2PO3 groups feature a dou- 
bly bonded oxygen in addition to the hydroxyl group. In 



these instances, binding is also facilitated by simultane- 
ous interaction of the methane with the 0-H bond and 
the lone-pair site, as seen in Fig. [3ljd)&(f). Examination 
of the charge density for these groups reveals that binding 
is accompanied by charge transfer from the lone pairs of 
the doubly bonded oxygen to the 0-H bond of the other 
oxygen, enhancing the polarity of the 0-H bond. The 
resulting configuration exploits the tetrahedral geometry 
of the adsorbate to maximize the specific interaction with 
the functional group. 



IV. CARBON DIOXIDE ADSORPTION 

We also calculated binding of CO2 on the functional- 
ized ZGNRs. DFT-f-LAP was used for these values, since 
it gave adequate results for methane binding distance, en- 
ergy, geometry, and trend across functional groups with 
respect to MP2 and available experimental data. The 
results for ii^ads and the final relaxed binding geometries 
are shown in Table HIl Fig. [Sj and Fig. \6\ According to 
Table HIl binding of CO2 usually exceeds that of CH4 
by about a factor of two, both for the ZGNRs and for 
clean graphene. This is due primarily to the presence of 
a quadrupole moment in CO2. The energetic difference 
between IP and OP adsorption is larger for CO2, with 
OP adsorption again preferred. Nevertheless, the order- 
ing of i?ads across the functional groups for both OP and 
IP adsorption is very similar for both gases. 

A comparison of the adsorption strength of CO2 versus 
CH4 for a given ZGNR offers insight into the relative se- 
lectivity of the substrate (last column of Table HI)) . Here, 
we find that the unfunctionalized (hydrogen-passivated) 
ZGNR shows an exceptionally high difference between 
CO2 and CH4 binding strengths: for DFT+LAP, OP 
(IP) binding of CO2 demonstrates a 257% (162%) in- 
crease over CH4, compared to the graphene baseline in- 
crease of 108%. This suggests that significantly enhanced 
selectivity for CO2 over CH4 may be achieved in graphitic 
systems with high concentrations of hydrogenated ex- 
posed edges. This is applicable to natural gas purifica- 
tion, in which efficient, selective CO2 removal from CH4 
streams is an important step. Nevertheless, edge con- 
centrations of oxygen- and nitrogen-containing functional 
groups should be kept sufficiently low, since these groups 
often decrease selectivity. Moreover, one should keep 
in mind that thermodynamic competition will favor sur- 
face sites over the fully hydrogen-passivated zigzag edge, 
meaning the realization of any practical benefit would 
likely require very high edge-carbon concentrations, as 
well as high gas pressures. 

For IP binding (Fig. [S|), the specific CO2 orientations 
vary from group to group. The OH, COOH, CH3, and 
NO2 ZGNRs favor alignment of the molecular C2 axis 
with the functional group (AT); the NH2 and unfunc- 
tionalized ribbon have the Coo axis aligned with the func- 
tional group (AP). In addition, most of the groups have 
the Coo axis of the IP CO2 molecule lying parallel to 




J^ 



(b)OH 



(c) NHg 





(d) COOH 



FIG. 5. Final adsorption geometries for IP adsorption of CO2 
on (a) an unfunctionalized ZGNR; as well as ZGNRs func- 
tionalized with (b) OH, (c) NH2, (d) COOH, (e) CH3, and 
(f) NO2. Configurations were obtained using the DFT-I-LAP 
method. For the OH- and COOH-functionalized ribbons, the 
CO2 molecule is rotated about its C2 axis by ~ 45° with 
respect to the ZGNR plane. Additional interactions of the 
CO2 molecule with nearby passivation hydrogens are evident 
in (a), (b), (c), and (f). Color scheme: H (cyan), C (yellow), 
O (red), N (green). The GO2 carbon is shown in gray. 






(a) H [none] 



(b)OH 



(c) NH2 





(e) CH3 




(f) NO2 



(g) graphene 



FIG. 6. Final adsorption geometries for OP adsorption of 
CO2 on (a) an unfunctionalized ZGNR; as well as ZGNRs 
functionalized with (b) OH, (c) NH2, (d) COOH, (e) CH3, 
and (f) NO2; and (g) on clean graphene. Configurations were 
obtained using the DFT-I-LAP method. Color scheme: H 
(cyan), C (yellow), O (red), N (green). The GO2 carbon is 
shown in gray. 
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TABLE II. Final relaxed geometries and adsorption energies 
iSads of CO2 adsorbed on functionalized and unfunctional- 
ized ZGNRs. The final orientation and position of the CO2 
molecule is given alongside the adsorption energy in parenthe- 
ses (see text for description of acronyms). For IP adsorption, 
dmin represents the distance from the CO2 carbon atom to 
the nearest functional-group atom, whose identity is given in 
parentheses. For OP adsorption, d^ is the distance of the ad- 
sorbate from the ribbon basal plane. For each value of -Eads, 
the last column shows the percent change with respect to the 
equivalent result for methane adsorption. 



Group 


Site 


Quantity 


DFT+LAP 


vs. CH4 


None (H) 


IP 


Sads (kJ/mol 


11.8 (AT) 
3.61 (H) 


-hl62% 




OP 


-Eads (kJ/mol 
d, (A) 


18.9 (AP,E) 
3.20 


+257% 


OH 


IP 


Sads (kJ/mol 

Ctmin (A) 


14.0 (AP) 
2.96 (H) 


-f57% 




OP 


-Eads (kJ/mol 
dz (A) 


20.8 (AP,E) 
3.09 


N/A 


NH2 


IP 


Eads (kJ/mol 

"min (A) 


17.8 (AT) 
3.34 (H) 


-F82% 




OP 


-Eads (kJ/mol 
d, (A) 


25.1 (AP,E) 
3.17 


+96% 


COOH 


IP 


-Bads (kJ/mol 

rfmin (A) 


19.9 (AP) 
2.78 (H) 


+101% 




OP 


-Eads (kJ/mol 
dz (A) 


21.2 (AP,E) 
3.23 


+86% 


CH3 


IP 


Sads (kJ/mol 

Ctmin (A) 


6.0 (AP) 
2.98 (H) 


+ 114% 




OP 


EJads (kJ/mol 
d, (A) 


22.3 (AP,E) 
3.21 


+ 125% 


NO2 


IP 


Sads (kJ/mol 


17.0 (AP) 
3.03 (0) 


+87% 




OP 


-Eads (kJ/mol 

d. (A) 


20.7 (AP,E) 
3.30 


+82% 


Graphene 


OP 


EJads (kJ/mol 

d. (A) 


23.1 (AP) 
3.14 


+108% 



the ZGNR plane. Exceptions are the OH- and COOH- 
functionalized ZGNRs, for which the Coo axis forms a 
nearly 45° angle with the plane of the ribbon. 

As we saw for the methane case, the ZGNRs function- 
alized with smaller polar functional groups (OH, NH2, 
NO2) feature significant secondary adsorbate interactions 
with neighboring edge hydrogens for IP binding. This is 
particularly true for NH2, where the AP binding geome- 
try results in the CO2 molecule being centered far along 
the ribbon axis (Fig.ISjc)). Also, the center carbon in the 
CO2 molecule aligns approximately along the polar bond 
in the NH2, COOH, and NO2 ZGNRs (Fig.[Si;c),(d),(f)), 
again similar to methane. Neighboring edge hydrogens 
are much more accessible on ZGNRs than on function- 



alized benzene, due to the latter's ring geometry; ac- 
cordingly, our geometries for ZGNRs with smaller polar 
groups differ from those reported in Ref. |6^ for CO2 on 
functionalized benzene. 

When looking at OP adsorption on the ZGNRs 
(Fig. ini, the location of the CO2 carbon, as well as its 
axial orientation, tends to mimic the graphene configura- 
tion. On graphene, the CO2 molecule sits parallel to the 
plane, such that the oxygens lie above the centers of two 
consecutive rings with the carbon above the center of a 
C-C bond (Fig. IH^g)). Similarly, for OP binding on the 
functionalized ZGNRs, the CO2 molecule always lies fiat 
in the ZGNR plane, with the CO2 carbon sitting above 
the ZGNR edges and never above the functional group 
atoms. This is even true for the unfunctionalized ZGNR 
(Fig. [DJa)). In addition, the molecule is always centered 
between two carbons on a zigzag ribbon edge, with one of 
its oxygens sitting above a ZGNR carbon ring, much like 
the graphene case. These observations echo our general 
conclusions for OP methane binding, where the preferred 
orientation on graphene significantly influenced the final 
methane configuration. The location of the remaining 
CO2 oxygen that does not sit above a carbon ring center 
varies. For NH2, CH3, and the unfunctionalized ZGNR 
(Fig. inKa),(c),(e)), it shows significant interaction with 
the functional group (or H, in the case of the unfunction- 
alized ribbon), resulting in alignment of the Coo axis with 
the primary axis of the group. For other ribbons, direct 
interaction is generally avoided. This choice of configu- 
ration may have some correlation with binding strength, 
since NH2 and CH3 are also the strongest-binding OP 
groups (see Table iHl). Note that unlike the methane case, 
DFT+LAP was able to identify a stable OP binding site 
for CO2 on the OH-functionafized ZGNR. 



V. DISCUSSION 

For both IP and OP adsorption, chemical functional- 
ization of the edges can enhance gas binding with respect 
to the unfunctionalized hydrogen-terminated ZGNR. 
This can be seen in Tables U and [III as well as in the 
summarized results for both gases using DFT+LAP, pre- 
sented in Fig. [T] The calculated degree of enhancement 
depends on the specific method and functionalization. 
Nevertheless, increases are universally demonstrated for 
all functionalizations except IP binding of CH3 when 
comparing within a given geometry (IP or OP). This sug- 
gests that exposed edges in carbon substrates, which are 
nominally rather inert to gas adsorption when passivated 
by hydrogen, can become binding sites when properly 
functionalized. 

For results in Tables U and [Tl] other than uncorrected 
DFT-LDA, OP binding is clearly preferred over IP bind- 
ing. For some groups, the differences between the two 
binding sites are quite large: CH3, for instance, demon- 
strates relatively strong OP binding in the DFT+LAP 
and vdW-DF descriptions but is by far the weakest- 



^ 20- 

-a 
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NO, 




-CO, 



-CH4 



\ ' \ ' 1 ' r 

12 3 4 

Dipole Moment Magnitude (Debye) 



FIG. 7. Calculated DFT+LAP adsorption energy £ads for 
CH4 (squares) and CO2 (circles) on all ZGNRs as a function 
of the permanent dipole moment of the correspondingly func- 
tionalized benzene substrate. Values for IP (solid line, open 
symbols) and OP (dashed line, filled symbols) adsorption are 
shown. The dotted lines track the average of the IP and OP 
values. Arrows on the right-hand side indicate the calculated 
values for OP adsorption on a clean graphene sheet. 



binding IP group. For others, such as CO OH, the bind- 
ing geometry has only a minor effect on i?ads in the van 
der Waals-corrected DFT schemes. However, it is worth 
emphasizing that even though OP binding seems to be 
energetically preferred, both geometries are relevant for 
the understanding and design of nanostructured carbon 
materials. This is because certain microstructures may 
preclude one binding geometry or the other. For instance, 
it is reasonable to assume that graphitic terraces may fa- 
vor IP binding, so as to permit simultaneous interaction 
with the terrace-attached functional group as well as the 
underlying surface layer. Similarly, both IP and OP bind- 
ing might be exhibited by a single adsorbate at a bilayer 
step edge, where it could interact simultaneously with the 
chemical functional group at the edge of one layer via IP 
binding, and at the edge of the second layer via OP bind- 
ing. In the case of interactions with functional groups in 
crystalline nanomaterials such as MOFs or COFs, steric 
factors could also favor either IP or OP adsorptioui^ 

It should be noted that the observed binding enhance- 
ments for the polar functionalizations (especially NH2, 
NO2, and COOH) over the hydrogen-passivated ZGNR 
apply only to edge binding. It may be argued that a 
more apt comparison for application purposes would be 
to a clean planar graphene surface, against which the 
associated increases in -Eads are effectively negligible for 
results other than uncorrected DFT-LDA (see Fig. [3 for 
example). However, the fact that edge adsorption can 
be made competitive with surface adsorption also means 
that edge sites become viable binding sites, thereby en- 
hancing the effective active surface area of the substrate. 
In nanoporous materials where edges are expected to 
comprise a significant fraction of available carbon binding 
sites (e.g., aerogels'^'^^), this should translate to enhanced 



capacity for gas uptake. Perhaps more significantly, a 
naive comparison of edge and surface binding belies the 
importance of the microstructure of the pore walls in 
nanoporous carbon, where edges and surfaces can be ex- 
pected to coexist. In these instances, one may assume 
the effect on the gas binding will be additive, with the 
functionalized edges providing additional stabilization at 
the pore wall beyond existent surface-site binding. 

The similarity in the results for CH4 and CO2 , despite 
differences in specific binding orientations and molecular 
symmetries, points to the importance of the substrate 
rather than the adsorbate in determining the binding 
trends. It is reasonable to assume that similar trends 
could be found for other physisorbed gas species. Note 
that the agreement between the CH4 and CO2 trends is 
particularly strong for IP binding, where the direct in- 
teraction with the functional group is likely to dominate 
relative to contributions from the bulk of the ribbon. 

In order to understand the specific origin of the ob- 
served trends in i?ads in the case IP binding, we have 
carried out a Lowdin population (LP) analysis on the 
DFT-I-LAP results to extract partial charges on the 
atoms. Changes in these quantities can be tracked to 
estimate the charge redistribution induced by adsorp- 
tion, both for the functional group and for the adsorbate. 
We point out that charges derived from an LP analysis 
are meaningful only in a relative sense; we use them as 
an approximate indicator of correlation between induced 
charges and observed adsorption trends. The LP charges 
induced upon IP adsorption of CH4 on the functionalized 
ZGNRs (OH, NH2, COOH, CH3, and NO2) are listed 
in Tables IIIII and IIVI for the functional group and the 
methane molecule, respectively. As expected for weak- 
binding systems, the magnitudes of the induced charges 
are small. The LP analysis shows that the charge distri- 
bution of the methane distorts slightly upon adsorption, 
resulting in induced methane multipole moments that 
interact electrostatically with the functionalized ZGNR 
(see Table lIVp . The magnitudes of the charges show 
that methane charge asymmetry is significantly larger 
for COOH, OH, NH2, and NO2 than for CH3, as ex- 
pected based on the values of i?ads in Table HI Also, in 
each of the stronger-binding cases, one of the methane 
C-H bonds is significantly more polarized than the oth- 
ers. Close examination reveals that for the OH, NH2, 
and COOH groups, this always involves the hydrogen 
atom that is farthest from the functional group, whereas 
for the strongly electron-withdrawing NO2 group, it in- 
volves the hydrogen that is closest to the group. The in- 
duced LP charges on the functional group atoms, shown 
in Table IIIII show similar results, with larger changes 
occurring for COOH, NH2, and OH than for the weak- 
binding CH3 . The calculations in Table IIIII also demon- 
strate that the largest changes are found for atoms within 
the polar 0-H and N-H bonds of the COOH, OH, and 
NH2 groups, further confirming that these polar bonds 
are primary contributors to methane binding. However, 
the NO2 group is a notable exception to these trends. 
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TABLE III. Induced charges (units of lO^^e) on the 
functional-group atoms after IP methane adsorption on a 
ZGNR, based on LP analysis. Results are obtained using 
the DFT+LAP method. (01, HI) refers to pairs of atoms in 
0-H bonds^ 

OH NH2 COOH CH3 NO2 



C -3.0 +2.7 

01 +4.5 +10.8 -1.3 

02 (=0) -3.5 -1.0 

N -4.5 -4.5 

HI -8.9 -8.4 -8.1 -3.3 
H2 -13.9 -2.8 

H3 -2.8 



functionalized benzene molecule possesses no dipole mo- 
ment, any exhibited polarization arises solely from the 
presence of the functional groups. Although OP and IP 
binding show somewhat different behavior, the average of 
-Eads for the two geometries correlates with the dipole mo- 
ment up to ^^2 Debye. -Eads saturates beyond this point, 
with the highly polar NO2 group showing no additional 
enhancement. Note that there is no obvious dependence 
of E'ads on whether the attached group is electron with- 
drawing or electron donating in nature. In other words, 
the physical amount of charge present in the it system is 
not as important for determining adsorption strength as 
is its arrangement, manifest in the induced polarization. 



TABLE IV. Induced charges (units of 10~ e) on atoms within 
methane after IP adsorption on a functionalized ZGNR, based 
on LP analysis. Results are obtained using the DFT+LAP 
method. 

OH NH2 COOH CH3 NO2 

C -9.9 -7.5 -18.1 +2.0 -1.9 

HI +13.4 +10.4 +12.6 +1.7 +3.0 

H2 +3.7 -3.1 +7.7 -3.6 +1.4 

H3 -6.8 -3.5 -4.8 -4.2 -7.1 

H4 -6.3 -5.6 -5.0 -4.5 -14.4 



exhibiting very weak induced LP charges on the func- 
tional group despite a rather large ii^ads and significant 
induced methane polarization. This could simply be an 
artifact of the LP partitioning scheme — for instance, no 
LP change would be observed if the charge redistribu- 
tion remains largely localized to the oxygen lone pairs. 
For this reason, it may be desirable to use a partitioning 
method that has a rigorous connection to the induced 
polarization, such as a methodology based on maximally 
localized Wannier functions. ^^ 

Although the LP analysis permits estimation of 
methane dipole moments, it is more difficult to ex- 
tract a measure of the polarization induced on the ex- 
tended ZGNRs upon adsorption. This is because of 
technical issues that arise when computing dipole mo- 
ments of extended structures, given that electronic po- 
larization along the repeating directions is ill defined 
within periodic boundary conditions!^ Alternatively, if 
we can assume that the dominant interactions between 
the methane and the substrate are local, we can use func- 
tionalized benzene molecules to estimate the dipole mo- 
ments induced on the functionalized ZGNRs. Since the 
molecular systems are isolated, there are no complica- 
tions in computing dipole moments. In Fig. [71 we plot 
the DFT+LAP E^ds for IP and OP gas binding on func- 
tionalized ZGNRs as a function of the magnitude of the 
permanent dipole moment of the corresponding benzene 
systems in the absence of an adsorbate, as calculated 
from DFT+LAP benzene models. Note that since an un- 



VL CONCLUSIONS 

In summary, we have shown that the binding of 
methane to graphitic carbon edges can be significantly 
modified by the presence of chemical functional groups. 
The degree of modification is found to be a function not 
only of the chemical group but also of the gas adsorp- 
tion geometry. Polar groups, including COOH, NH2, 
NO2, and H2PO3 are found to be promising candidates 
for enhancing CO2 and CH4 capacity by activating ex- 
posed edges to introduce additional binding sites. This 
has important implications for improving uptake in car- 
bon structures that are known to have high quantities of 
edge carbons, including nanoporous activated complexes, 
carbon aerogels, and crystalline covalent framework ma- 
terials. 

In contrast, by far the highest selectivity for CO2 
over CH4 was found for the unfunctionalized, hydrogen- 
passivated zigzag edge. In this particular case, the CO2 
binding preference was significantly higher than for bulk 
graphene. This implies that higher edge concentrations 
could be beneficial for natural gas purification using 
nanoporous carbon, provided the substrates contain rel- 
atively low concentrations of other functional groups. 
Nevertheless, since the hydrogen-passivated zigzag edge 
was also among the weakest adsorbers, combining high 
selectivity with high uptake for maximum separation ef- 
fectiveness remains a design challenge. 

The observed trends in binding energies with func- 
tional group identity are largely preserved across the 
DFT-LDA, DFT+LAP, and vdW-DF formulations of 
the DFT exchange-correlation functional if one considers 
only in-plane or only out-of-plane binding sites. Never- 
theless, comparisons between in-plane and out-of-plane 
binding energies are not qualitatively consistent across 
the methods, suggesting that one should use caution in 
directly comparing results based on very different binding 
geometries. 

We expect that our work will stimulate experimental 
efforts to confirm our results, by processing activated car- 
bons to favor the presence of these chemical functional 
groups and looking for a concomitant change in methane 
uptake. This should constitute a realistic step towards 
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